﻿using System;
using System.Collections.Generic;
using SpiceNet.Circuits;
using SpiceNet.Diagnostics;
using MathNet.Numerics.LinearAlgebra;
using System.Diagnostics;
using SpiceNet.Models;

namespace SpiceNet.Simulations
{
    /// <summary>
    /// An abstract class definition for an analysis
    /// </summary>
    public abstract class Analysis
    {
        /// <summary>
        /// An event that is called when new data is available
        /// </summary>
        public event ExportDataEventHandler ExportData;

        /// <summary>
        /// Types of domains
        /// </summary>
        public enum DomainTypes
        {
            NoDomain = 0x00,
            TimeDomain = 0x01,
            FrequencyDomain = 0x02,
            SweepDomain = 0x03
        }

        /// <summary>
        /// The name of the analysis
        /// </summary>
        public string Name { get; protected set; }

        /// <summary>
        /// The description of the analysis
        /// </summary>
        public string Description { get; protected set; }

        /// <summary>
        /// Gets a list of jobs for this type of analysis
        /// </summary>
        public IEnumerable<Job> Jobs => jobs;
        private List<Job> jobs;

        /// <summary>
        /// Gets or sets the domain of the analysis
        /// </summary>
        public DomainTypes Domain { get; protected set; }

        /// <summary>
        /// Do IC or not?
        /// </summary>
        public bool DoIc { get; protected set; }

        /// <summary>
        /// Initialize the analysis
        /// </summary>
        public virtual void Init(Circuit ckt)
        {
            return;
        }

        /// <summary>
        /// Run the function of the analysis
        /// </summary>
        public abstract void Func(Circuit ckt, bool reset);

        /// <summary>
        /// Export data
        /// </summary>
        /// <param name="ckt"></param>
        protected void Dump(Circuit ckt)
        {
            ExportDataArgs args = new ExportDataArgs(ckt);

            // Export data if any event is attached
            ExportData?.Invoke(this, args);
        }

        /// <summary>
        /// Private variables
        /// </summary>
        private Type jobtype = null;

        /// <summary>
        /// Constructor
        /// </summary>
        /// <param name="jobtype"></param>
        protected Analysis(Type jobtype)
        {
            this.jobtype = jobtype;
            jobs = new List<Job>();
        }

        /// <summary>
        /// Add a job to the analysis
        /// </summary>
        /// <param name="job"></param>
        public void AddJob(Job job)
        {
            if (job.GetType() != jobtype)
                throw new InvalidJobTypeException(this);
            jobs.Add(job);
        }

        /// <summary>
        /// Remove a job
        /// </summary>
        /// <param name="job"></param>
        public void RemoveJob(Job job) => jobs.Remove(job);

        /// <summary>
        /// Find the operating point of the circuit
        /// </summary>
        /// <param name="ckt"></param>
        /// <param name="firstmode"></param>
        /// <param name="continuemode"></param>
        /// <param name="iterlim"></param>
        public void Op(
            Circuit ckt,
            Circuit.Modes firstmode,
            Circuit.Modes continuemode,
            int iterlim)
        {

            ckt.Mode = firstmode;
            if (!ckt.Config.NoOpIter)
            {
                if (Iterate(ckt, iterlim))
                    return;
            }

            // No convergence
            // try Gmin stepping
            if (ckt.Config.NumGminSteps > 1)
            {
                ckt.Mode = firstmode;
                CircuitWarning.Warning(this, "Starting Gmin stepping");
                ckt.DiagGmin = ckt.Config.Gmin;
                for (int i = 0; i < ckt.Config.NumGminSteps; i++)
                    ckt.DiagGmin *= 10.0;
                for (int i = 0; i <= ckt.Config.NumGminSteps; i++)
                {
                    ckt.IsCon = true;
                    if (!Iterate(ckt, iterlim))
                    {
                        ckt.DiagGmin = 0.0;
                        CircuitWarning.Warning(this, "Gmin step failed");
                        break;
                    }
                    ckt.DiagGmin /= 10.0;
                    ckt.Mode = continuemode;
                }
                ckt.DiagGmin = 0.0;
                if (Iterate(ckt, iterlim))
                    return;
            }

            // No we'll try source stepping
            if (ckt.Config.NumSrcSteps > 1)
            {
                ckt.Mode = firstmode;
                CircuitWarning.Warning(this, "Starting source stepping");
                for (int i = 0; i <= ckt.Config.NumSrcSteps; i++)
                {
                    ckt.SrcFact = ((double)i) / (double)ckt.Config.NumSrcSteps;
                    if (!Iterate(ckt, iterlim))
                    {
                        ckt.SrcFact = 1.0;
                        ckt.CurrentAnalysis = AnalysisType.DoingTran;
                        CircuitWarning.Warning(this, "Source stepping failed");
                        return;
                    }
                }
                ckt.SrcFact = 1.0;
                return;
            }

            // Failed
            throw new CircuitException("Could not determine operating point");
        }

        /// <summary>
        /// Iterate using a circuit Load function
        /// NOTE: This function can often not converge (eg. Transient analysis with too big timesteps)
        /// Be careful though when porting from Spice though, it uses 0 when convergence succeeded! This 
        /// function returns TRUE if it converges. To make matters worse, Spice usually calls the return 
        /// value "converged", which gives the illusion that it is != 0 in case of convergence.
        /// </summary>
        /// <param name="ckt">The circuit</param>
        /// <param name="iterlim">The maximum amount of iterations to be tried</param>
        /// <returns>True if the iteration succeeded</returns>
        protected bool Iterate(Circuit ckt, int iterlim)
        {
            int iterno = 0;
            bool ipass = false;

            // Ignore operating condition point, just use the solution as-is
            if (ckt.Mode.HasFlag(Circuit.Modes.TranOp) && ckt.Mode.HasFlag(Circuit.Modes.Uic))
            {
                // Swap solution and old solution
                ckt.State.SwapSolutions();
                Load(ckt);
                return true;
            }

            if (!ckt.State.Initialized)
                ckt.State.Init(ckt);

            // Perform iteration
            while (true)
            {
                ckt.IsCon = true;

                if (!ckt.Mode.HasFlag(Circuit.Modes.InitPred))
                {
                    try
                    {
                        Load(ckt);
                        iterno++;
                    }
                    catch (CircuitException)
                    {
                        iterno++;
                        ckt.Stat.NumIter = iterno;
                        throw;
                    }

                    // Solve the equation (thank you Math.NET)
                    ckt.Stat.SolveTime.Start();
                    ckt.State.Solve();
                    ckt.Stat.SolveTime.Stop();

                    // Exceeded maximum number of iterations
                    if (iterno > iterlim)
                    {
                        ckt.Stat.NumIter += iterno;
                        return false;
                    }

                    if (ckt.IsCon && iterno != 1)
                        ckt.IsCon = IsConvergent(ckt);
                    else
                        ckt.IsCon = false;
                }

                if ((ckt.Mode & Circuit.Modes.InitFloat) != 0)
                {
                    if ((ckt.Mode & Circuit.Modes.Dc) != 0 && ckt.HadNodeset)
                    {
                        if (ipass)
                            ckt.IsCon = false;
                        ipass = false;
                    }

                    // We have convergence, return
                    if (ckt.IsCon)
                    {
                        ckt.Stat.NumIter += iterno;
                        return true;
                    }
                }
                else if ((ckt.Mode & Circuit.Modes.InitJct) != 0)
                    ckt.Mode = (ckt.Mode & ~Circuit.Modes.InitF) | Circuit.Modes.InitFix;
                else if ((ckt.Mode & Circuit.Modes.InitFix) != 0)
                {
                    if (ckt.IsCon)
                        ckt.Mode = (ckt.Mode & ~Circuit.Modes.InitF) | Circuit.Modes.InitFloat;
                    ipass = true;
                }
                else if ((ckt.Mode & Circuit.Modes.InitSmSig) != 0)
                    ckt.Mode = (ckt.Mode & ~Circuit.Modes.InitF) | Circuit.Modes.InitFloat;
                else if ((ckt.Mode & Circuit.Modes.InitTran) != 0)
                    ckt.Mode = (ckt.Mode & ~Circuit.Modes.InitF) | Circuit.Modes.InitFloat;
                else if ((ckt.Mode & Circuit.Modes.InitPred) != 0)
                    ckt.Mode = (ckt.Mode & ~Circuit.Modes.InitF) | Circuit.Modes.InitFloat;
                else
                {
                    ckt.Stat.NumIter += iterno;
                    throw new CircuitException("Could not find flag");
                }

                // We need to do another iteration, swap solutions with the old solution
                ckt.State.SwapSolutions();
            }
        }

        /// <summary>
        /// Load the circuit in the Rhs vector and Matrix
        /// </summary>
        protected void Load(Circuit ckt)
        {
            // Start the stopwatch
            ckt.Stat.LoadTime.Start();

            // Clear rhs and matrix
            ckt.State.ClearLoad();

            // Load all devices
            foreach (Type t in ckt.ModelTypes)
            {
                Model[] models = ckt[t];
                for (int i = 0; i < models.Length; i++)
                    models[i].Load(ckt);
            }

            // Check modes
            if ((ckt.Mode & Circuit.Modes.Dc) != 0)
            {
                // Consider doing nodeset & ic assignments
                if ((ckt.Mode & (Circuit.Modes.InitJct | Circuit.Modes.InitFix)) != 0)
                {
                    // Do nodesets
                    for (int i = 0; i < ckt.Nodes.Count; i++)
                    {
                        CircuitNode node = ckt.Nodes[i];
                        if (node.Nodeset.Given)
                        {
                            if (ZeroNoncurRow(ckt.State.Matrix, ckt.Nodes, node.Number))
                            {
                                ckt.State.Rhs[node.Number] = 1.0e10 * node.Nodeset;
                                ckt.State.Matrix[node.Number, node.Number] = 1.0e10;
                            }
                            else
                            {
                                ckt.State.Rhs[node.Number] = node.Nodeset;
                                ckt.State.Solution[node.Number] = node.Nodeset;
                                ckt.State.Matrix[node.Number, node.Number] = 1.0;
                            }
                        }
                    }
                }
                if ((ckt.Mode & Circuit.Modes.TranOp) != 0 && (ckt.Mode & Circuit.Modes.Uic) == 0)
                {
                    for (int i = 0; i < ckt.Nodes.Count; i++)
                    {
                        CircuitNode node = ckt.Nodes[i];
                        if (node.IC.Given)
                        {
                            if (ZeroNoncurRow(ckt.State.Matrix, ckt.Nodes, node.Number))
                            {
                                ckt.State.Rhs[node.Number] = 1.0e10 * node.IC;
                                ckt.State.Matrix[node.Number, node.Number] = 1.0e10;
                            }
                            else
                            {
                                ckt.State.Rhs[node.Number] = node.IC;
                                ckt.State.Solution[node.Number] = node.IC;
                                ckt.State.Matrix[node.Number, node.Number] = 1.0;
                            }
                        }
                    }
                }
            }

            // Keep statistics
            ckt.Stat.LoadTime.Stop();
        }

        /// <summary>
        /// Test convergence
        /// </summary>
        /// <param name="ckt"></param>
        private bool IsConvergent(Circuit ckt)
        {
            // Check convergence for each node
            for (int i = 0; i < ckt.Nodes.Count; i++)
            {
                CircuitNode node = ckt.Nodes[i];
                double n = ckt.State.Solution[node.Number];
                double o = ckt.State.OldSolution[node.Number];
                if (node.Type == NodeType.Voltage)
                {
                    double tol = ckt.Config.RelTol * Math.Max(Math.Abs(n), Math.Abs(o)) + ckt.Config.VoltTol;
                    if (Math.Abs(n - o) > tol)
                    {
                        return false;
                    }
                }
                else
                {
                    double tol = ckt.Config.RelTol * Math.Max(Math.Abs(n), Math.Abs(o)) + ckt.Config.AbsTol;
                    if (Math.Abs(n - o) > tol)
                    {
                        // Convergence failed
                        // ckt.TroubleNode = i;
                        return false;
                    }
                }
            }

            // Check for convergence with the models
            foreach (Type t in ckt.ModelTypes)
            {
                Model[] models = ckt[t];
                for (int i = 0; i < models.Length; i++)
                {
                    // Returns true if error
                    if (!models[i].IsConvergent(ckt))
                        return false;
                }
            }

            // Convergence succeeded
            return true;
        }

        /// <summary>
        /// Reset the row to 0.0 and return true if the row is a current equation
        /// </summary>
        /// <param name="matrix"></param>
        /// <param name="nodes"></param>
        /// <param name="rownum"></param>
        /// <returns></returns>
        protected bool ZeroNoncurRow(Matrix<double> matrix, NodeList nodes, int rownum)
        {
            bool currents = false;
            for (int n = 0; n < nodes.Count; n++)
            {
                CircuitNode node = nodes[n];
                double x = matrix[rownum, node.Number];
                if (x != 0.0)
                {
                    if (node.Type == NodeType.Current)
                        currents = true;
                    else
                        matrix[rownum, node.Number] = 0.0;
                }
            }
            return currents;
        }
    }

    /// <summary>
    /// A class for exporting data
    /// </summary>
    public class ExportDataArgs : EventArgs
    {
        /// <summary>
        /// Get the circuit
        /// </summary>
        public Circuit Ckt { get; }

        /// <summary>
        /// Constructor
        /// </summary>
        /// <param name="ckt"></param>
        public ExportDataArgs(Circuit ckt) : base()
        {
            Ckt = ckt;
        }
    }

    /// <summary>
    /// An event handler for exporting data for the current analysis
    /// </summary>
    /// <param name="a"></param>
    /// <param name="ckt"></param>
    public delegate void ExportDataEventHandler(object sender, ExportDataArgs e);

    /// <summary>
    /// Types of analysis
    /// </summary>
    public enum AnalysisType
    {
        DoingNone = 0x00,
        DoingDcOp = 0x01,
        DoingTrCv = 0x02,
        DoingAc = 0x04,
        DoingTran = 0x08
    }
}
